home *** CD-ROM | disk | FTP | other *** search
/ CU Amiga Super CD-ROM 19 / CU Amiga Magazine's Super CD-ROM 19 (1998)(EMAP Images)(GB)[!][issue 1998-02].iso / CUCD / Programming / LEDA / source / src / plane / _sweep_segments.c < prev    next >
Encoding:
C/C++ Source or Header  |  1994-11-16  |  13.4 KB  |  587 lines
  1. /*******************************************************************************
  2. +
  3. +  LEDA  3.1c
  4. +
  5. +
  6. +  _sweep_segments.c
  7. +
  8. +
  9. +  Copyright (c) 1994  by  Max-Planck-Institut fuer Informatik
  10. +  Im Stadtwald, 6600 Saarbruecken, FRG     
  11. +  All rights reserved.
  12. + 
  13. *******************************************************************************/
  14.  
  15.  
  16. //------------------------------------------------------------------------------
  17. //
  18. // SWEEP_SEGMENTS
  19. //
  20. // a plane sweep algorithm for computing line segment intersections
  21. //
  22. // SWEEP_SEGMENTS(list<segment>& L1, list<segment>& L2, GRAPH<point,int> G);
  23. //
  24. // computes the planar subdivision created by the (directed) straight
  25. // line segments of L1 and L2
  26. // 
  27. // G = (V,E) with
  28. //
  29. // 1. each node v in V contains a point of intersectoin between two segments
  30. //    (from L1 or L2)
  31. //
  32. // 2. each edge e = (v,w) corresponds to a
  33. //
  34. //
  35. //------------------------------------------------------------------------------
  36.  
  37.  
  38. #include <LEDA/sweep_segments.h>
  39. #include <LEDA/sortseq.h>
  40. #include <LEDA/prio.h>
  41. #include <math.h>
  42.  
  43.  
  44. #define EPS  0.00001
  45. #define EPS2 0.0000000001
  46.  
  47. class S_POINT;
  48. class S_SEGMENT;
  49. typedef S_POINT* S_point;
  50. typedef S_SEGMENT* S_segment;
  51.  
  52. enum S_point_type {Cross=0,Rightend=1,Leftend=2}; 
  53.  
  54. class S_POINT
  55. {
  56.   friend class S_SEGMENT;
  57.  
  58.   S_segment seg;
  59.   int     kind;
  60.   double  x;
  61.   double  y;
  62.  
  63.   public:
  64.  
  65.   S_POINT(double a,double b)  
  66.   { 
  67.     x=a; y=b; seg=0; kind=Cross;
  68.    }
  69.  
  70.  
  71.  
  72.   S_POINT(point p)         
  73.   { 
  74.     x=p.xcoord();y=p.ycoord();seg=0;kind=Cross;
  75.    }
  76.  
  77.  
  78.   LEDA_MEMORY(S_POINT);
  79.  
  80.   friend double    get_x(S_point);
  81.   friend double    get_y(S_point);
  82.   friend int       get_kind(S_point);
  83.   friend S_segment get_seg(S_point);
  84.  
  85.   friend bool intersection(S_segment, S_segment, S_point&);
  86. };
  87.  
  88. inline double    get_x(S_point p)    { return p->x; }
  89. inline double    get_y(S_point p)    { return p->y; }
  90. inline int       get_kind(S_point p) { return p->kind; }
  91. inline S_segment get_seg(S_point p)  { return p->seg; }   
  92.  
  93.  
  94.  
  95. static int compare(const S_point& p1, const S_point& p2)
  96. { if (p1==p2) return 0;
  97.  
  98.   double diffx = get_x(p1) - get_x(p2);
  99.   if (diffx >  EPS2 ) return  1;
  100.   if (diffx < -EPS2 ) return -1;
  101.  
  102.   int    diffk = get_kind(p1)-get_kind(p2);
  103.   if (diffk != 0) return diffk;
  104.  
  105.   double diffy = get_y(p1) - get_y(p2);
  106.   if (diffy >  EPS2 ) return  1;
  107.   if (diffy < -EPS2 ) return -1;
  108.  
  109.   return 0;
  110.  
  111. }
  112.  
  113.  
  114. class S_SEGMENT
  115. {
  116.   S_point startpoint;
  117.   S_point endpoint;
  118.   double  slope;
  119.   double  yshift;
  120.   node  left_node;
  121.   int   orient;
  122.   int   color;
  123.   int   name;
  124.  
  125.  
  126.   public:
  127.  
  128.   S_SEGMENT(S_point, S_point,int,int);     
  129.  
  130.  ~S_SEGMENT() { delete startpoint; delete endpoint; }     
  131.  
  132.   LEDA_MEMORY(S_SEGMENT);
  133.  
  134.   friend S_point get_startpoint(S_segment);
  135.   friend S_point get_endpoint(S_segment);
  136.   friend double get_slope(S_segment);
  137.   friend double get_yshift(S_segment);
  138.   friend int  get_orient(S_segment);
  139.   friend int  get_color(S_segment);
  140.   friend int  get_name(S_segment);
  141.   friend node get_left_node(S_segment);
  142.   friend void set_left_node(S_segment, node);
  143.  
  144.   friend bool intersection(S_segment, S_segment, S_point&);
  145. };
  146.  
  147.  
  148. inline S_point get_startpoint(S_segment seg)     { return seg->startpoint; }
  149. inline S_point get_endpoint(S_segment seg)       { return seg->endpoint; }
  150. inline double get_slope(S_segment seg)           { return seg->slope; }
  151. inline double get_yshift(S_segment seg)          { return seg->yshift; }
  152. inline int  get_orient(S_segment seg)            { return seg->orient; }
  153. inline int  get_color(S_segment seg)             { return seg->color; }
  154. inline int  get_name(S_segment seg)              { return seg->name; }
  155. inline node get_left_node(S_segment seg)         { return seg->left_node; }
  156. inline void set_left_node(S_segment seg, node v) { seg->left_node = v; }
  157.  
  158.  
  159.  
  160. S_SEGMENT::S_SEGMENT(S_point p1,S_point p2,int c, int n)    
  161.   {
  162.     left_node  = nil;
  163.     color      = c;
  164.     name       = n;
  165.  
  166.     if (compare(p1,p2) < 0)
  167.      { startpoint = p1; 
  168.        endpoint = p2; 
  169.        orient = 0;
  170.       }
  171.     else
  172.      { startpoint = p2; 
  173.        endpoint = p1; 
  174.        orient = 1;
  175.       }
  176.  
  177.     startpoint->kind = Leftend; 
  178.     endpoint->kind = Rightend; 
  179.     startpoint->seg = this; 
  180.     endpoint->seg = this;
  181.  
  182.     if (endpoint->x != startpoint->x)
  183.     {
  184.       slope = (endpoint->y - startpoint->y)/(endpoint->x - startpoint->x);
  185.       yshift = startpoint->y - slope * startpoint->x;
  186.  
  187.       startpoint->x -= EPS;
  188.       startpoint->y -= EPS * slope;
  189.       endpoint->x += EPS;
  190.       endpoint->y += EPS * slope;
  191.     }
  192.     else //vertical segment
  193.     { startpoint->y -= EPS;
  194.       endpoint->y   += EPS;
  195.      }
  196.   }
  197.  
  198.  
  199. static double x_sweep;
  200. static double y_sweep;
  201.  
  202.  
  203. static int compare(const S_segment& s1, const S_segment& s2)
  204. {
  205.   double y1 = get_slope(s1)*x_sweep+get_yshift(s1);
  206.   double y2 = get_slope(s2)*x_sweep+get_yshift(s2);
  207.  
  208.   double diff = y1-y2;
  209.   if (diff >  EPS2 ) return  1;
  210.   if (diff < -EPS2 ) return -1;
  211.  
  212.   if (get_slope(s1) == get_slope(s2)) 
  213.         return compare(get_x(get_startpoint(s1)), get_x(get_startpoint(s2)));
  214.  
  215.   if (y1 <= y_sweep+EPS2)
  216.         return compare(get_slope(s1),get_slope(s2));
  217.   else
  218.         return compare(get_slope(s2),get_slope(s1));
  219.  
  220. }
  221.  
  222.  
  223. void Print(S_segment& x) 
  224. { S_point s = get_startpoint(x);
  225.   S_point e = get_endpoint(x);
  226.   cout << 
  227.     string("[(%f,%f)----(%f,%f)]",get_x(s),get_y(s), get_x(e),get_y(e));
  228. }
  229.  
  230.  
  231. static priority_queue<seq_item,S_point> X_structure;
  232.  
  233. static sortseq<S_segment,pq_item> Y_structure;
  234.  
  235.  
  236. bool intersection(S_segment seg1,S_segment seg2, S_point& inter)
  237. {
  238.   if (seg1->slope == seg2->slope)
  239.     return false;
  240.   else
  241.   { 
  242.     double cx = (seg2->yshift - seg1->yshift) / (seg1->slope - seg2->slope);
  243.  
  244.     if (cx <= x_sweep) return false;
  245.  
  246.     if (seg1->startpoint->x > cx || seg2->startpoint->x > cx ||
  247.         seg1->endpoint->x < cx || seg2->endpoint->x < cx ) return false;
  248.  
  249.     inter = new S_POINT(cx,seg1->slope * cx + seg1->yshift);
  250.  
  251.     return true;
  252.   }
  253. }
  254.  
  255.  
  256. static pq_item Xinsert(seq_item i, S_point p) 
  257. { 
  258.   return X_structure.insert(i,p);
  259. }
  260.  
  261. static S_point Xdelete(pq_item i) 
  262. {
  263.   S_point p = X_structure.inf(i);
  264.   X_structure.del_item(i);
  265.   return p;
  266. }
  267.  
  268.  
  269. static node New_Node(GRAPH<point,int>& G,double x, double y )
  270. { return G.new_node(point(x,y)); }
  271.  
  272. static void New_Edge(GRAPH<point,int>& G,node v, node w, S_segment l )
  273. { if (get_orient(l)==0)
  274.        G.new_edge(v,w,get_color(l));
  275.   else G.new_edge(w,v,get_color(l));
  276. }
  277.  
  278.  
  279. void handle_vertical_segment(GRAPH<point,int>& SUB, S_segment l)
  280. { 
  281.   S_point p = new S_POINT(get_x(get_startpoint(l)),get_y(get_startpoint(l)));
  282.   S_point q = new S_POINT(get_x(get_endpoint(l)),get_y(get_endpoint(l)));
  283.  
  284.   S_point r = new S_POINT(get_x(p)+1,get_y(p));
  285.   S_point s = new S_POINT(get_x(q)+1,get_y(q));
  286.  
  287.   S_segment bot = new S_SEGMENT(p,r,0,0);
  288.   S_segment top = new S_SEGMENT(q,s,0,0);
  289.  
  290.   seq_item bot_it = Y_structure.insert(bot,nil);
  291.   seq_item top_it = Y_structure.insert(top,nil);
  292.   seq_item sit;
  293.  
  294.   node u,v,w;
  295.   S_segment seg;
  296.   
  297.  
  298.   for(sit=Y_structure.succ(bot_it); sit != top_it; sit=Y_structure.succ(sit))
  299.   { seg = Y_structure.key(sit);
  300.  
  301.     double cross_y = (get_slope(seg) * get_x(p) + get_yshift(seg));
  302.  
  303.     v = get_left_node(seg);
  304.     if (v==nil)
  305.     { w = New_Node(SUB,get_x(p),cross_y);
  306.       set_left_node(seg,w);
  307.      }
  308.     else
  309.     { double vx = SUB[v].xcoord();
  310.       if ( vx < get_x(p)-EPS2) 
  311.       { w = New_Node(SUB,get_x(p),cross_y);
  312.         New_Edge(SUB,v,w,seg);
  313.         set_left_node(seg,w);
  314.        }
  315.       else w = v;
  316.      }
  317.  
  318.     u = get_left_node(l);
  319.     if (u!=nil && u!=w) New_Edge(SUB,u,w,l);
  320.     set_left_node(l,w);
  321.  
  322.    }
  323.   
  324.   delete l;
  325.   delete top;
  326.   delete bot;
  327.     
  328.   Y_structure.del_item(bot_it);
  329.   Y_structure.del_item(top_it);
  330.  
  331.  }
  332.  
  333.  
  334. void SWEEP_SEGMENTS(const list<segment>& S, list<point>& P)
  335. { GRAPH<point,int> G;
  336.   list<segment> S0;
  337.  
  338.   SWEEP_SEGMENTS(S,S0,G);
  339.  
  340.   node v;
  341.   forall_nodes(v,G) P.append(G[v]);
  342.  
  343. }
  344.  
  345. void SWEEP_SEGMENTS(const list<segment>& L1, 
  346.                     const list<segment>& L2, GRAPH<point,int>& SUB)
  347. {
  348.   S_point    p,inter;
  349.   S_segment  seg, l,lsit,lpred,lsucc,lpredpred;
  350.   pq_item  pqit,pxmin;
  351.   seq_item sitmin,sit,sitpred,sitsucc,sitpredpred;
  352.  
  353.  
  354.   int count=1;
  355.  
  356.   //initialization of the X-structure
  357.  
  358.   segment s;
  359.  
  360.   forall(s,L1) 
  361.    { S_point p = new S_POINT(s.start());
  362.      S_point q = new S_POINT(s.end());
  363.      seg = new S_SEGMENT(p,q,0,count++);
  364.      Xinsert(nil,get_startpoint(seg));
  365.    }
  366.  
  367.   count = -1;
  368.  
  369.   forall(s,L2) 
  370.    { S_point p = new S_POINT(s.start());
  371.      S_point q = new S_POINT(s.end());
  372.      seg = new S_SEGMENT(p,q,1,count--);
  373.      Xinsert(nil,get_startpoint(seg));
  374.    }
  375.  
  376.  
  377.   count = 0;
  378.  
  379.   x_sweep = -MAXINT;
  380.   y_sweep = -MAXINT;
  381.  
  382.  
  383.   while(!X_structure.empty())
  384.   {
  385.     pxmin = X_structure.find_min();
  386.     p = X_structure.inf(pxmin);
  387.  
  388.     sitmin = X_structure.key(pxmin);
  389.  
  390.     Xdelete(pxmin);
  391.  
  392.  
  393.  
  394.     if (sitmin == nil) //left endpoint
  395.     { 
  396.  
  397.       l = get_seg(p); 
  398.  
  399.       x_sweep = get_x(p);
  400.       y_sweep = get_y(p);
  401.  
  402.  
  403.       if (get_x(p) == get_x(get_endpoint(l)))
  404.         handle_vertical_segment(SUB,l);
  405.       else
  406.       {
  407.  
  408. /*
  409.       sit = Y_structure.lookup(l);
  410.       if (sit!=nil)  
  411.            error_handler(1,"plane sweep: sorry, overlapping segments");
  412. */
  413.  
  414.  
  415.       sit = Y_structure.insert(l,nil);
  416.  
  417.       Xinsert(sit,get_endpoint(l));
  418.  
  419.       sitpred = Y_structure.pred(sit);
  420.       sitsucc = Y_structure.succ(sit);
  421.  
  422.       if (sitpred != nil) 
  423.       { if ((pqit = Y_structure.inf(sitpred)) != nil)
  424.           delete Xdelete(pqit);
  425.  
  426.         lpred = Y_structure.key(sitpred);
  427.  
  428.         Y_structure.change_inf(sitpred,nil);
  429.  
  430.         if (intersection(lpred,l,inter))
  431.             Y_structure.change_inf(sitpred,Xinsert(sitpred,inter));
  432.       }
  433.  
  434.  
  435.       if (sitsucc != nil)
  436.       { lsucc = Y_structure.key(sitsucc);
  437.         if (intersection(lsucc,l,inter))
  438.            Y_structure.change_inf(sit,Xinsert(sit,inter));
  439.       }
  440.      } /* else if vertical */
  441.  
  442.     }
  443.     else if (get_kind(p) == Rightend)
  444.          //right endpoint
  445.          { 
  446.            x_sweep = get_x(p);
  447.            y_sweep = get_y(p);
  448.  
  449.            sit = sitmin;
  450.  
  451.            sitpred = Y_structure.pred(sit);
  452.            sitsucc = Y_structure.succ(sit);
  453.  
  454.            S_segment seg = Y_structure.key(sit);
  455.  
  456.            Y_structure.del_item(sit);
  457.  
  458.            delete seg;
  459.  
  460.            if((sitpred != nil)&&(sitsucc != nil))
  461.            {
  462.              lpred = Y_structure.key(sitpred);
  463.              lsucc = Y_structure.key(sitsucc);
  464.              if (intersection(lsucc,lpred,inter))
  465.                 Y_structure.change_inf(sitpred,Xinsert(sitpred,inter));
  466.            }
  467.          }
  468.          else // point of intersection
  469.          { 
  470.            node w = New_Node(SUB,get_x(p),get_y(p));
  471.  
  472.            count++;
  473.  
  474.            /* Let L = list of all lines intersecting in p 
  475.  
  476.               we compute sit     = L.head();
  477.               and        sitpred = L.tail();
  478.  
  479.               by scanning the Y_structure in both directions 
  480.               starting at sitmin;
  481.  
  482.            */
  483.  
  484.            /* search for sitpred upwards from sitmin: */
  485.  
  486.            Y_structure.change_inf(sitmin,nil);
  487.  
  488.            sitpred = Y_structure.succ(sitmin);
  489.  
  490.  
  491.            while ((pqit=Y_structure.inf(sitpred)) != nil)
  492.            { S_point q = X_structure.inf(pqit);
  493.              if (compare(p,q) != 0) break; 
  494.              X_structure.del_item(pqit);
  495.              Y_structure.change_inf(sitpred,nil);
  496.              sitpred = Y_structure.succ(sitpred);
  497.             }
  498.  
  499.  
  500.            /* search for sit downwards from sitmin: */
  501.  
  502.            sit = sitmin;
  503.  
  504.            seq_item sit1;
  505.            
  506.            while ((sit1=Y_structure.pred(sit)) != nil)
  507.            { pqit = Y_structure.inf(sit1);
  508.              if (pqit == nil) break;
  509.              S_point q = X_structure.inf(pqit);
  510.              if (compare(p,q) != 0) break; 
  511.              X_structure.del_item(pqit);
  512.              Y_structure.change_inf(sit1,nil);
  513.              sit = sit1;
  514.             }
  515.  
  516.  
  517.  
  518.            // insert edges to p for all S_segments in sit, ..., sitpred into SUB
  519.            // and set left node to w 
  520.  
  521.            lsit = Y_structure.key(sitpred);
  522.  
  523.            node v = get_left_node(lsit);
  524.            if (v!=nil) New_Edge(SUB,v,w,lsit);
  525.            set_left_node(lsit,w);
  526.  
  527.            for(sit1=sit; sit1!=sitpred; sit1 = Y_structure.succ(sit1))
  528.            { lsit = Y_structure.key(sit1);
  529.  
  530.              v = get_left_node(lsit);
  531.              if (v!=nil) New_Edge(SUB,v,w,lsit);
  532.              set_left_node(lsit,w);
  533.             }
  534.  
  535.            lsit = Y_structure.key(sit);
  536.            lpred=Y_structure.key(sitpred);
  537.            sitpredpred = Y_structure.pred(sit);
  538.            sitsucc=Y_structure.succ(sitpred);
  539.  
  540.  
  541.            if (sitpredpred != nil)
  542.             { 
  543.               lpredpred=Y_structure.key(sitpredpred);
  544.  
  545.               if ((pqit = Y_structure.inf(sitpredpred)) != nil)
  546.                 delete Xdelete(pqit);
  547.  
  548.               Y_structure.change_inf(sitpredpred,nil);
  549.  
  550.  
  551.               if (intersection(lpred,lpredpred,inter))
  552.                 Y_structure.change_inf(sitpredpred,
  553.                                        Xinsert(sitpredpred,inter));
  554.              }
  555.  
  556.  
  557.            if (sitsucc != nil)
  558.             {
  559.               lsucc=Y_structure.key(sitsucc);
  560.  
  561.               if ((pqit = Y_structure.inf(sitpred)) != nil)
  562.                 delete Xdelete(pqit);
  563.                  
  564.               Y_structure.change_inf(sitpred,nil);
  565.  
  566.               if (intersection(lsucc,lsit,inter))
  567.                   Y_structure.change_inf(sit,Xinsert(sit,inter));
  568.              }
  569.  
  570.  
  571. // reverse the subsequence sit, ... ,sitpred  in the Y_structure
  572.  
  573.            x_sweep = get_x(p);
  574.            y_sweep = get_y(p);
  575.  
  576.            Y_structure.reverse_items(sit,sitpred);
  577.  
  578.           delete p;
  579.  
  580.          } // intersection
  581.  
  582.   }
  583.  
  584.   X_structure.clear();
  585.  
  586. }
  587.